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Abstract: Direct geo-referencing is an efficient methodology for the fast acquisition of 
3D spatial data. It requires the fusion of spatial data acquisition sensors with navigation 
sensors, such as Global Navigation Satellite System (GNSS) receivers. In this contribution, 
we consider an integrated GNSS navigation system to provide estimates of the position 
and attitude (orientation) of a 3D laser scanner. The proposed multi-sensor system (MSS) 
consists of multiple GNSS antennas rigidly mounted on the frame of a rotating laser scanner 
and a reference GNSS station with known coordinates. Precise GNSS navigation requires 
the resolution of the carrier phase ambiguities. The proposed method uses the multivariate 
constrained integer least-squares (MC-LAMBDA) method for the estimation of rotating 
frame ambiguities and attitude angles. MC-LAMBDA makes use of the known antenna 
geometry to strengthen the underlying attitude model and, hence, to enhance the reliability 
of rotating frame ambiguity resolution and attitude determination. The reliable estimation 
of rotating frame ambiguities is consequently utilized to enhance the relative positioning of 
the rotating frame with respect to the reference station. This integrated (array-aided) method 
improves ambiguity resolution, as well as positioning accuracy between the rotating frame 
and the reference station. Numerical analyses of GNSS data from a real-data campaign 
confirm the improved performance of the proposed method over the existing method. In 
particular, the integrated method yields reliable ambiguity resolution and reduces position 



Sensors 2014, 14 



12716 



standard deviation by a factor of about 0.8, matching the theoretical gain of a/3/4 for two 
antennas on the rotating frame and a single antenna at the reference station. 

Keywords: global navigation satellite system (GNSS); attitude determination; multivariate 
constrained integer least- squares (MC-LAMBDA); carrier phase ambiguity resolution; direct 
geo-referencing; laser scanner 



1. Introduction 

The acquisition and interpretation of three-dimensional (3D) spatial data are important assets for 
scientific and industrial applications, such as 3D city modeling, facility management, construction 
engineering, navigation and forensic investigations. Direct geo-referencing, which does not require 
dedicated control points, is an efficient methodology for the fast acquisition of 3D spatial data by means 
of a 3D laser scanner. It can be performed either using additional backsight targets [1-4] or using 
external sensors [5-8]. The latter requires the fusion of spatial data acquisition sensors and navigation 
sensors, such as Global Navigation Satellite System (GNSS) sensors. In this contribution, we consider 
an integrated GNSS navigation system to provide estimates of the position and attitude (orientation) of 
a 3D laser scanner. 

The use of GNSS for geo-referencing has been explored in various studies. Direct geo-referencing is 
demonstrated using GNSS integrated with inertial sensors [9,10] and a digital compass [5]. In this work, 
we explore pure a GNSS-based navigation solution as in [8,1 1]. In [11], a single rotating antenna is used 
to provide a post-processing navigation solution. As in [7,8], the proposed multi-sensor system (MSS) 
consists of multiple GNSS antennas rigidly and symmetrically mounted on the frame of a rotating laser 
scanner and a reference GNSS station with known coordinates. 

The proposed method uses the multivariate constrained integer least-squares (MC-LAMBDA) 
method [12-17] for the estimation of rotating frame ambiguities and attitude angles. MC-LAMBDA 
makes use of known antenna geometry to strengthen the underlying attitude model, enabling reliable 
instantaneous ambiguity resolution and attitude determination of the rotating frame. The reliable 
estimation of rotating frame ambiguities is consequently utilized to enhance the positioning of the 
rotating frame. This array-aided positioning method [15,18-20] naturally yields the estimates of the 
rotating frame center (centroid of antennas' reference points) and improves ambiguity resolution, as well 
as the positioning accuracy of the relative position between the rotating frame and the reference station. 

The numerical studies considered in this contribution include performance analyses of the proposed 
method with GNSS data from two real data campaigns. Comparison studies using epoch-by-epoch 
processing and filtering confirm the improved performance of the proposed method over the existing 
method from [8]. This contribution is organized as follows: Section 2 describes the multi-sensor 
system considered and defines the problem. Section 3 describes our attitude determination and filtering 
approaches for the rotating frame. Section 4 describes the array-aided positioning and filtering methods 
for the positioning of the rotating frame. Section 5 presents real data analyses demonstrating the 
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improved performance of the proposed method. Finally, Section 6 contains the summary and conclusions 
of this contribution. 

2. Background 

The multi-sensor system (MSS) considered for geo-referencing in this contribution consists of a laser 
scanner and two GNSS antennas/receivers. As shown in Figure la, the laser scanner is the core sensor 
of the MSS, which rotates about its vertical axis with a constant angular velocity. GNSS receivers are 
connected to two eccentric GNSS antennas, which are mounted such that the centroid of the antenna 
reference points (ARPs) coincides with the scanner rotating axis. In addition to these GNSS receivers, it 
is assumed to have a nearby reference GNSS station with a known position (Figure lb). During the data 
acquisition, the MSS makes a complete 360 degree rotation about its vertical axis, collecting both laser 
and GNSS measurements, which are synchronized through a GNSS receiver event marker. 

Figure 1. Geo-referencing scenario, (a) The MSS formed by a laser scanner (blue) and two 
eccentrically-mounted GNSS antennas (green) [8]; and (b) the reference GNSS station with 
a known position. 




(a) (b) 

The objective of the navigation system is to provide the position (the centroid of ARPs) and 
the pointing direction (heading) of the laser scanner. In [8], standard real-time kinematic (RTK) 
positioning [21] is used to estimate individual rotating antenna positions, and then, a constrained 
nonlinear filtering method, in particular an extended Kalman filter, is used to obtain the above 
parameters. In this contribution, we use constrained integer least-squares (Section 3) and array-aided 
positioning (Section 4), enabling improved ambiguity resolution and improved positioning accuracy. 
In the following sections, we formulate a more general problem, estimating attitude angles and relative 
position between two platforms with multiple GNSS antennas/receivers, which enables us to demonstrate 
the potential of array-aided positioning. As shown in the following sections, array-aided positioning 
utilizes the reliable estimation of rotating frame ambiguities, which are obtained from array processing 
(MC-LAMBDA), improving the estimation of the relative position of the rotating frame with respect to 
the given reference station. 
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3. Attitude Determination 

This section describes the platform processing involving attitude determination for a small-sized array 
of GNSS receivers/antennas with a known local body frame antenna geometry. First, the multi-baseline 
attitude model is introduced using the multivariate formulation of [12]. This formulation makes frequent 
use of the Kronecker product and the vec-operator [22]. Then, we include the local body frame 
antenna-geometry and show how the constrained attitude model can be solved in a step-wise manner. 

3.1. The Multivariate Model 

Let us consider the A;-th platform equipped with a set of n k + 1 antennas simultaneously tracking 
m + 1 satellites on / frequencies. The set of linearized double difference (DD) GNSS phase and code 
observations obtained on the n k baselines formed by these antennas at an observation epoch forms a 
multivariate Gauss-Markov model [12,19]: 

E{Y k ) = AZ k + GB k , Z k e Z fmxnk (1) 
D(vec(Y k )) = Q Ykyk = P nk ® Q yy , B k G M 3xrifc (2) 

where E(-) and D(-) denote the expectation and dispersion operator, <g> denotes the Kronecker product, 
Y k = [y k , . . . , y k k ] the 2fm x n k matrix of n k linearized (observed-minus-computed) DD observation 
vectors y k , Z k = [z k , . . . , z k k ] the fm x n k matrix of n k unknown DD integer ambiguity vectors z k , 
B k = [b k , ... ,b k } the 3 x n k matrix of n k unknown baseline vectors bj, G the 2fm x 3 geometry 
matrix that contains the unit line-of-sight vectors, A the 2fm x fm matrix that links the DD data to the 
integer ambiguities and P nk and Q yy the known matrices of order n k x n k and 2fm x 2fm, respectively. 
Here, vec(-) denotes the vec-operator, which transforms a matrix into a vector by stacking the columns 
of the matrix, one underneath the other. Note that, for the simplicity of the formulation, we assumed 
that all receivers/antennas track the same set of satellites. However, this restriction is relaxed in the 
software implemented using MATLAB. Since the unit line-of-sight vectors of two antennas to the same 
satellite on a short baseline considered in this work (<10 km) are the same for all practical purposes, the 
geometry matrix G is considered the same for different platforms, as well as for the between-platform 
baseline at a given time instant. 

For the stochastic model, we assumed that all receivers/antennas have similar (noise) characteristics. 
However, the results in the following are also applicable for dissimilar receivers/antennas [19]. The 
correlation matrix P nk takes care of the correlation that follows from the fact that the n k baselines share 
the observations from the reference receiver. For similar receivers/antennas, it is given as: 

Pn k = ^ + e n k ^ nk ) (3) 

with I nk the identity matrix of size n k and e nfc the n k -vector of ones. Matrix Q yy takes care of the 
precision of the phase and code data and is given as: 

Q yy = blockdiag (Q l , . . . , Q f ) (4) 

where: 



Q f = 2 x blockdiag (Qf :p , 



(5) 
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with Q f:p = D T m Q' f:p D m , Q ft9 = DlQ' f .^D m , Q' f:t = diag [(a}. t ) 2 , . . . , (a^ +1 ) 2 J , D T m = [-e m I m ] 
the single difference operator, "blockdiag" referring to the block diagonal matrix formed by given 
arguments and "diag" referring to the diagonal matrix formed by given arguments. The factor 
two in Equation (5) is due to the between-receiver difference of similar receivers. We assume 
elevation-dependent noise characteristics [23] for the undifferenced observables. That is, the standard 
deviation of the undifferenced observable can be written as: 

4* = a /:*o ( X + a f--t 0 ex P ) (6) 

where 9 s is the elevation angle of satellite s and crf-t 0 , o-f:t 0 and 9f :tQ are the elevation-dependent 
model parameters. 

3.2. The Body -Frame Antenna-Geometry as Multivariate Constraints 

The strength of the above model can be improved by including information about the geometry of the 
antenna configuration. The known body-frame antenna-geometry can be included into the above model 
through the parametrization: 

B k = R k B k , R k e ® 3xqk (7) 

with the unknown 3 x q k orthogonal matrix R k (R kT R k = I qk ), 0 3X9fc denoting the set of orthogonal 
matrices of size 3 x ^ and the known q k x n k matrix Bq = [b k 1; . . . , &o n J describing the known 
geometry of the antenna configuration in the body frame. Here, q k is the degree of geometrical 
independence of the GNSS baselines, for example, q k = 1 for co-linearly installed antennas, q k = 2 
for co-planarly installed antennas and q k = 3 for antennas not installed in a single plane. For q k = 3, 
R k is related to the Euler attitude angles $ =[4>6 as follows: 



R($) 



-s e s^ce c^pce 



(8) 



with 0 the heading, 9 the elevation, ip the bank and where s a = sin(a) and c a = cos(a). Note that for 
q < 3, only the first q columns of R are defined. For example, for a linear antenna array (q = 1), only 
the first column is defined, and hence, only heading and elevation are estimable. For q > 1 (an array 
with more than two antennas that are not in a straight line), all three angles are estimable. 

Substitution of Equation (7) into Equation (1) leads to the constrained GNSS attitude model [19,24]: 

E(Y k ) = AZ k + GR k B k Z k e Z fmxnk (9) 
D(vec(F fc )) = Q YkYk =P nk ®Qyy _R fe G 0 3x<7fe (10) 

Our objective is to solve for the attitude matrix R k in a least-squares sense, thereby taking the integer 
constraint on matrix Z k E Z^ mxnfc and the orthonormality constraint on matrix R k e 0 3X9fe into account. 
Hence, the least-squares minimization problem that will be solved reads: 

min \\vec(Y k - AZ k -GR k B k )\\ 2 n (11) 
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with || • ||q = (-) T Q This is a mixed integer nonlinear least-squares problem that does not permit 
a closed-form solution. We now describe how Equation (11) can be solved. 

3.3. The Real-Valued Float Solution 

The float solution is defined as the solution of Equation (11) without the constraints. When we ignore 
the integer constraint on Z k and the orthonormality constraint on R k , the float solutions Z k and R k and 
their variance-covariance matrices are obtained from solving the system of normal equations: 



with: 



Qz k Z k Qz k R k 

Qk k z k Qfi k R k 



Qz k z k Qz k h k 

Q ' R k Z k Q ' R k R k 



-1 


vec(Z k ) 




vec(R k ) 



(12) 



(AlQ- Y l Yk &k) \ At 



In k ® A T 

B k ® G T 



The ^^-constrained solution of R k and its variance-covariance matrix can be obtained from the float 
solution as follows: 

vec(i? fc (Z fc )) = vec(R k ) - Q R k z k Q^z kVec ~ Z& ) 

Q R k (Z k )R k (Z k ) = Qk k R k ~ Q R k Z k Q £k £kQ Z k R k 

= (B k P ni T l B kT y l <g) {G T Q7yyGy 1 
Using the above estimators, the original problem in Equation (11) can be decomposed as: 

i2 



(13) 



(14) 



min 1 1 vec (Y k - AZ k - GR k B k ) \ 

z k ez /m x n k ,_r* eo 3 x q k 



vec ( E 



+ min 

,k v k Z k £Zf mxn k 



vec Z k -Z k 



Q) - L. - I. 



+ min 

Rk eQ 3xq k 



vec[R k (Z ky 



R 1 



Q R k (Z k )R k {Z k ) 



(15) 



with E k = Y k — AZ k — GR k B-Q being the matrix of least-squares residuals. Note that the first term on 
the right-hand side is constant, as it does not depend on the unknown matrices Z k and R k . 

3.4. The Integer Ambiguity Solution 

Based on the orthogonal decomposition (15), the multivariate constrained integer minimization can 
be formulated as: 



Z k = arg min C k (Z k ) 

z k ezf mxn k 



(16) 



where: 



C k (Z k " 



vec(Z k - Z k ) 



+ 



vec (R k (Z k ) -R k (Z k ) 



Qr 



(17) 



R k (Z k )R k (Z k ) 
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with: 

(18) 



R k (Z k ) = arg min vec ( R k (Z k ) - R k 



'R k (Z k )R k (Z k ) 



The ambiguity objective function C k (Z k ) is the sum of two coupled terms: the first weighs the 
distance from the float ambiguity matrix Z k to the nearest integer matrix Z k in the metric of Q z k z k, while 
the second weighs the distance from the conditional float solution R k (Z k ) to the nearest orthonormal 
matrix R k in the metric of Q R k (z k )R k (z k r Unlike with the standard LAMBDA method [25], the search 
space of the above integer minimization problem is non-ellipsoidal, due to the presence of the second 
term in C k (Z k ). This second term is a consequence of having the orthonormality constraints rigorously 
included. The evaluation of C k (Z k ) requires the computation of a nonlinear constrained least-squares 
problem (18) for every integer matrix in the search space. In the MC-LAMBDA method, this problem is 
mitigated through the use of easy-to-evaluate bounding functions [17]. 



3.5. The Ambiguity Resolved Attitude Solution 

Finally, we obtain the integer ambiguity-resolved attitude solution by substituting Z k into 
Equation (13), thus giving R k (Z k ). The sought-for attitude angles i3 k (Z k ) are then given by 
the reparametrized solution of Equation (18). Using a first order approximation, the formal 
variance-co variance matrix of the attitude angles is given by: 

where J R kjk is the Jacobian of -d k (R k ). 



3.6. Attitude Filtering 

The epoch-by-epoch MC-LAMBDA attitude solution is further processed using an unscented Kalman 
filter (UKF) [26]. For a leveled platform (i.e., for small 9 and ip), the kinematic equations of the attitude 
angles are given as [27]: 



a>i = + vf_ x 



(20) 



where the state vector a. 



<f>i fa 9i 9i ipi ipi 



state transition matrix F is given as: 



consists of attitude angles and angular rates, and the 



F = h 



1 T 
0 1 



(21) 



where T is the sampling interval. The process noise vf_ x has a zero mean normal distribution with 
variance-covariance matrix Q v <*v a ,i-i, which is given as: 



Q 



v a v a ,i-l 



diag ([<rj, erg, crj]) ® 



T 3 /3 T 2 /2 
T 2 /2 T 



(22) 
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with ct0, 00 and 0^ the process noise standard deviations. The observation model reads: 

Q = h(on) + wf 



(23) 



with Q given by (ir^Z* )J at epoch i. The nonlinear observational function h(oii) is defined by 
Equation (8), and the observation noise wf is assumed to have a zero mean normal distribution with 



covariance matrix Q w a wa i , which is given by Q n 



$ k (Z k )R k (Z k ) 



at epoch i. 



The use of the above constant-velocity model [28] reflects the fact that the frame is rotating at a 
constant rate. For the two-antenna scenario considered in real-data analyses (Figure 2), only heading and 
elevation angles are estimable. Hence, a reduced state space model consisting of only heading, elevation 
and their rates is used in Section 5. The recursive filter is initialized with two-point initialization [28] 

3 3 

and propagated with process noise standard deviations of 0^ = 0.01°s~2 and a e = 0°s"2 (i.e., dead 
reckoning filtering for elevation constraining to horizontal ID rotation). 

Figure 2. Multi- sensor experiment set-up on the roof of the building of the Geodetic Institute 
(Messdach), Leibniz Universitat Hannover, Germany. The MSS is mounted on Pillar 5, 
while the reference station is located on Pillar 8. (a) Location; and (b) the MSS used in 
the experiments. 




(a) 



(b) 



4. Integrated Positioning 

This section describes the between-platform processing involving relative positioning between two 
platforms equipped with arrays of GNSS receivers/antennas. The array-aided positioning described in 
the following is a novel positioning concept improving between-platform positioning using an array 
of antennas on the platforms [15,18-20]. Unlike the formulations in [18-20], the formulation in this 
contribution explicitly considers different numbers of antennas on the reference and user platforms. 
Unlike the parameter space formulation in [15], the current contribution considers a simplified, 
double-difference observation space formulation elegantly demonstrating the advantages of array-aided 
positioning. First, the combined observation model for all independent baselines among all receivers 
on both platforms is described. Then, we describe attitude-bootstrapping, showing how platform arrays 
improve the between-platform baseline estimate. 
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Let us consider two platforms carrying rii + 1 and n 2 + 1 receivers/antennas. The functional and 
stochastic models for the between-platform baseline formed by the first antennas (pivot antennas) read: 



E(y 12 ) = Az 12 + Gb 12 z 12 eZ fm (24) 

D(Z/ 12 ) = Qyy (25) 

where y 12 is the between-platform double-difference observables, z 12 is the unknown between-platform 
double-difference ambiguities and b 12 is the unknown between-platform baseline. Note that atmosphere 
delays are not considered in this formulation, as troposphere delays and ionosphere delays can be ignored 
for the short baseline (<10 km) considered in this work. However, these atmosphere delays must be 
taken into account for general long baseline scenarios [19]. In standard positioning, the LAMBDA 
method yields the optimal estimates for the ambiguities and, hence, for the baseline. 



4.1. Integrated Between-Platform Model 

By combining between-platform observables in Equation (24) and platform array observables in 
Equation (9), the functional and stochastic models of the integrated system read: 



where y = [Y 1 Y 2 



E(y) = AZ + G^B 0 
D(vecC7")) = P®Q yy 
y 12 } is the combined observation matrix, ^ 



is the combined rotation matrices and between-platform baseline, Z 



[Z 1 



R 2 
Z 2 



z 



(26) 
(27) 

£12] G ]^3x( gi +g 2 +l) 
121 G %fmxn t is the 



combined ambiguity matrix, B 0 = blockdiag (B^, B%, 1) is the combined local geometry matrix and: 



(28) 



is the combined correlation matrix with n t = ni+n 2 + l. The above system consists of attitude models of 
both platforms with unknowns Z k and R k and the between-platform baseline model with unknowns z 12 
and b 12 . Even though these three subsystems do not have any parameter in common, they are correlated 
as in Equation (28), due to the use of common observations from pivot antennas. 





p 

1 rii 


0 ie 


p = 


0 


P --e 






--e T 1 




. 2°ni 





4.2. Attitude Bootstrapping 

The attitude bootstrapping method [18,19] uses a decorrelation technique to decouple the 
combined system in Equation (26), such that the subsystems still yield the optimal solution. Using 
decorrelation matrix: 



D 



The decorrelated system reads: 



0 



0 



1 712 



0 
0 



l p T p-1 I t p-1 I 
2 ni 1 ni 2 ri2 n 2 



< 2/m 



(29) 



EOT) 
D (vec(y')) 



AZ' + G^B 0 

P ® Qyy 



(30) 
(31) 
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where y' = [Y 1 Y 2 y' 12 ~\ is the decorrelated observation matrix, %^ = [R 1 R 2 b' 12 ] G 
^3x(gi+ g 2+i) [ s j-^g combined rotation matrices and between-platform baseline after decorrelation, 
ZJ = [Z 1 Z 2 z' 12 ] is the combined ambiguity matrix after decorrelation and: 



with: 



P' = blockdiag(P ni , P n2 , rj) (32) 



1 _ ni _ 1 "2 

y 12 = ^-^tE^ + ^ttE^ 2 (33) 

r=l r=l 



^^^ + ^7^^ (34) 

+ 1 ^ r n 2 + 1 ^ r 

r=l r=l 
1 ni n 2 



z '12 = .12 



ni 1 n 2 

'12 _ l12 x \^ nll.1 , 1 \^ D 2i2 



m + 1 ^ < n 2 + _ p=i 

77 2(m + l)(n 2 + l) 1 ^ 

This decorrelation keeps the platform processing intact as in Equation (16) and only alters the 
between-platform model. As a result of decorrelation, the ambiguities in Equation (34) may not be 
an integer. However, once platform ambiguities are determined reliably using MC-LAMBDA with 
decoupled platform models in Equation (30), the model for the between-platform baseline can be 
rearranged as: 

E(y" 12 ) = Az 12 + Gb' 12 z 12 EZ fm (37) 
D(/ 12 ) = VQyy (38) 



where: 



y" 12 = y vz + —rl^ Az r-—jl^ Az r (39) 

Tlx + 1 ^-f n 2 + 1 

r=l r=l 

Due to the reduction of variance-covariance by a factor of r\, this model yields improved ambiguity 
resolution and baseline estimation compared to the standard positioning model in Equation (24). That 
is, the use of array-aided positioning reduces the variance-covariance matrices of the float ambiguities 
and ambiguity-fixed baseline estimators by a factor of r\. For the rotating frame with two antennas and 
a single antenna at the reference station in Figure 1, the variance reduction factor is r\ = |. Note that 
the between-platform baseline estimate in Equation (35) corresponds to between-centroids of antenna 
arrays, naturally yielding the parameter of interest for the geo-referencing system in Figure 1. The 
unconstrained mixed-integer least-squares problem defined in Equations (37) and (38) can be solved 
efficiently using the LAMBDA method [25] providing ambiguity-fixed baseline estimate b' 12 (Z>') and 
associated variance-covariance matrix Q y 12 (z')b' 12 (z')- 



r=l r=l 
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4.3. Baseline Filtering 

The epoch-by-epoch baseline solution in Section 4.2 is further processed to obtain the filtered 
estimates for the center of the MSS (Figure la). Unlike the previous method in [8], which uses 
constrained nonlinear filtering for antenna positions, the integrated method in Section 4.2 yields 
estimates of the center position, which is assumed to be stationary for a rotating, leveled frame. Hence, 
dead reckoning (linear) Kalman filtering yields the filtered estimates for the stationary center position. 
The kinematic equation reads: 

Pi = Pi-i + vL (40) 

where the state vector = b'P consists of position components. The process noise has a zero mean 
normal distribution with variance-covariance matrix Q v p v p ti -x, which is given as: 

QvPvM-i = diag (VL Oby, ObS) (41) 

with db x , Ob y and Ob z the process noise standard deviations. Since the center position of the rotating 
frame is stationary, a dead reckoning filter is used (i.e., Ob x = (Jb y = &b z = 0, which is equivalent to the 
recursive least-squares estimation of constant parameter vector b 12 ). The observation model reads: 

6 = A + tuf (42) 

with £j given by b' 12 (Z') at epoch i. Observation noise wf is assumed to have a zero mean normal 
distribution with covariance matrix Q w p w Pj, which is given by Qb'^(z')b' 12 {z') at e P ocn i- 

5. Analyses 

For numerical analyses, we used the data from a static and a kinematic experiment on the roof of the 
Geodetic Institute (Messdach) building, Leibniz Universitat Hannover, Germany. The MSS is mounted 
on Pillar 5 (cf. Figure 2) and equipped with a terrestrial laser scanner (TLS) Z+F Imager 5006 and 
two LEIAX1202GG GNSS antennas about 0.6 m apart. These antennas are connected to two JAVAD 
TRE_G3TH DELTA GNSS receivers. The reference station is located on Pillar 8 (about 20 m from 
the MSS and equipped with a JAVAD TRE.G3TH DELTA GNSS receiver and a LEIAR25.R3 LEIT 
antenna. For the kinematic experiment, we also considered another reference station equipped with a 
LEICA GRX1200 GNSS receiver and a LEIAR25.R4 LEIT antenna and located about 6 km away from 
the MSS. 

The static experiment was conducted on 4 October 2011, for about six hours with the collection of 
GPS data at a rate of 1 Hz. In the kinematic experiment on 7 October 2011, we collected GPS data 
for five subsequent rotations (about an hour) at a rate of 10 Hz. For all of our analyses, we used the 
elevation-dependent stochastic model with the parameters given in Table 1 . 
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Table 1. Elevation-dependent stochastic model parameters defined in Equation (6). 







Code 




Phase 




Frequency (/) 


°7:po 


a f--Po 9f-Po 




a f:<po 






(cm) 


(deg) 


(mm) 




(deg) 


LI and L2 


15 


5 20 


1 


5 


20 



5.1. Static Data 

This section presents the analyses of the static data demonstrating the benefits of using the knowledge 
of antenna geometry for attitude determination and array-aided positioning. Figure 3 shows the 
satellite visibility: the skyplot, the number of satellites and the position dilution of precision (PDOP) 
values. With a 10° elevation cut-off, we have a moderate GPS satellite geometry (PDOP < 3) during 
the experiment. 

Figure 3. GPS satellite visibility during the static experiment with a 10° elevation cut-off. 
(a) The number of satellites and position dilution of precision (PDOP); and (b) skyplot. 



o 




Table 2 summarizes the empirical instantaneous integer ambiguity resolution success rate for attitude 
determination, indicating the advantage of using MC-LAMBDA for the case of single-frequency 
processing. Similarly, Table 3 demonstrates the improved success rate performance of the proposed 
array-aided positioning with two antennas/receivers on the frame. Note that further improvement is 
possible by having more antennas/receivers. 
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Table 2. Empirical instantaneous ambiguity resolution success rate (%) for attitude 
determination using static data (0.6-m baseline). 



Processing 


LAMBDA 


MC-LAMBDA 


Single-frequency 


90.4 


100 


Dual-frequency 


100 


100 



Table 3. Empirical instantaneous ambiguity resolution success rate (%) for baseline 
estimation using static data (20-m baseline). 



Processing 


Standard Positioning 


Array- Aided Positioning 


Single-frequency 


85.8 


90.0 


Dual-frequency 


100 


100 



Figure 4 depicts the scatter plot of the ambiguity-fixed attitude angles, while Figure 5 shows the 
plots (scatter plot of the horizontal components and time series of the down component) for the 
ambiguity-fixed baseline estimates. Tables 4 and 5 summarise the corresponding empirical standard 
deviations. Note that, once the ambiguities have been resolved, the precision of the attitude solution 
is driven by the high precision of the carrier phase observations [16]. That is, the accuracy of the 
unconstrained attitude solution (using the LAMBDA method) is comparable to that of the constrained 
solution (using the MC-LAMBDA method), provided that ambiguities are correctly fixed. However, the 
knowledge of the antenna geometry plays an important role by strengthening the underlying model and, 
hence, improving the ambiguity resolution (see Table 2). In the case of baseline estimation (Table 5), 
the proposed method yields slightly improved estimates, due to the integrated processing with an array 
of antennas. 

Figure 4. Scatter plot of the ambiguity-fixed attitude angles using epoch-by-epoch 
processing of the static data (0.6-m baseline). 
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Figure 5. North-east scatter plot of the ambiguity-fixed baseline estimation using static data 
with epoch-by-epoch processing (20-m baseline). 

15| ! 1 ! 1 ! 1 




Table 4. Empirical angular standard deviation (deg) for the attitude determination using the 
epoch-by-epoch processing of the static data (0.6-m baseline). 



Heading 


Elevation 


0.24 


0.38 



Table 5. Empirical position standard deviation (mm) for baseline estimation using the 
epoch-by-epoch processing of the static data (20-m baseline). 



Processing Method 


North 


East 


Up 


Standard positioning 


3.4 


2.1 


5.1 


Array-aided positioning 


2.8 


1.6 


4.4 



5.2. Kinematic Data 

This section presents the analyses of the kinematic data comparing the proposed method with the 
existing method. These dual-frequency GPS data analyses aim to compare the estimation of the 
parameters of interest for geo-referencing, namely the heading and center point of the rotating frame. 
Figure 6 shows the satellite visibility: the skyplot, the number of satellites and the PDOP values. With a 
10° elevation cut-off, we have a good GPS satellite geometry (PDOP ~ 2) during the experiment. 

Table 6 summarizes the root mean square (RMSE) values of the heading angle for all five rotations, 
where the ground truth is determined using the fact that the frame was rotating at a constant rate and 
synchronized though a GNSS receiver event marker. Since the precision of ambiguity-fixed attitude 
angles is driven by the high precision of the carrier phase observations, both the proposed method and 
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the previous method [8] have a similar angular accuracy. Filtering further improves the estimates. Based 
on these analyses, the achievable heading angular accuracy using a 0.6-m baseline on the rotating frame 
is about 0.2° RMSE. 

Figure 6. GPS satellite visibility during the kinematic experiment with a 10° elevation 
cut-off. (a) The number of satellites and PDOP; and (b) skyplot. 

o 
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Table 6. Heading root mean square (RMSE) (deg) for kinematic data (0.6-m baseline). 



Rotation 


Epoch-by-Epoch 


Filtering 


Previous 


Proposed Previous 


Proposed 




Method 


Method 


Method 


Method 


1 


0.24 


0.20 


0.20 


0.16 


2 


0.26 


0.23 


0.24 


0.20 


3 


0.23 


0.16 


0.21 


0.13 


4 


0.20 


0.15 


0.16 


0.11 


5 


0.20 


0.17 


0.16 


0.14 



Baseline estimation errors of epoch-by-epoch processing and filtering for the first rotation are depicted 
in Figure 7, while the 3D position RMSE values for all five rotations are reported in Table 7. Here, 
smoothing estimates based on all five rotations are considered as the ground truth. The apparent 
improved performance of the proposed method is due to the use of novel integrated processing. Since 
the proposed integrated method naturally yields the estimates of the center point (see Section 4.2), 
the proposed filtering has a simple and strong dynamic model compared to the nonlinear, constrained 
filtering of the existing method [8] and, hence, yields improved estimates. Based on these analyses, the 
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achievable position accuracy using a 20-m baseline with two antennas on the rotating frame is about 
3 mm RMSE. 

Figure 7. Estimation error (mm) for kinematic data (20-m baseline) using the proposed 
method: Rotation 1 . (a) North error; and (b) east error. 
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Finally, we considered the determination of the center point using a reference station at about 6 km 
away. The baseline estimation results for this long baseline using the proposed method are provided in 
Table 8. The significant increase of the position RMSE of this baseline compared to the short baseline 
case in Table 7 is due to the presence of residual atmosphere delays for this long baseline. Hence, the 
achievable position RMSE for this baseline is about 20 mm. 
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Table 7. 3D position root mean square (RMSE) (mm) for kinematic data (20-m baseline). 



Epoch-by-Epoch Filtering 

Rotation 

Previous Proposed Previous Proposed 

Method Method Method Method 



1 


6.1 


3.1 


3.5 


1.5 


2 


4.2 


2.9 


3.4 


1.6 


3 


5.5 


3.2 


2.9 


1.5 


4 


5.1 


3.4 


3.2 


1.5 


5 


5.3 


4.7 


3.4 


2.6 



Table 8. 3D position RMSE (mm) for kinematic data (6-km baseline) processing with the 
proposed method. 



Rotation 


Epoch-by-Epoch 


Filtering 


1 


12.4 


9.7 


2 


12.3 


10.4 


3 


15.3 


14.3 


4 


21.0 


20.0 


5 


13.1 


7.8 



6. Summary and Conclusions 

In this contribution, we explored the use of an array of GNSS antennas for attitude determination and 
relative positioning for direct geo-referencing. The MC-LAMBDA method exploits the known antenna 
geometry to improve the reliability of resolving rotating frame ambiguities and, hence, to improve the 
reliability of the rotating frame attitude determination. Furthermore, the reliable estimation of rotating 
frame ambiguities enables the strengthening of the estimation of the baseline between the rotating frame 
and a reference station. Our analysis includes epoch-by-epoch processing, as well as recursive filtering. 
We demonstrated the improved performance of the proposed method using data from two experiments 
with a prototype MSS representing a rotating frame. The use of constrained attitude determination 
and array-aided positioning increases the reliability (in terms of ambiguity resolution) and improves 
the achievable position accuracy. It enables instantaneous ambiguity resolution, which is immune to 
cycle slips, and, hence, enables instantaneous mapping. Furthermore, the reliability and accuracy can 
further be improved by employing more antennas on the rotating frame and at the reference station. 
With a sufficient number of low-cost GNSS receivers, the potential of instantaneous mobile mapping for 
low-cost applications will be explored in future studies. 
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